Breakdown of the Stokes-Einstein relation in two, three and four dimensions 
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The breakdown of the Stokes-Einstein (SE) relation between diffusivity and viscosity at low tem- 
peratures is considered to be one of the hallmarks of glassy dynamics in liquids. Theoretical analy- 
ses relate this breakdown with the presence of heterogeneous dynamics, and by extension, with the 
fragility of glass formers. We perform an investigation of the breakdown of the SE relation in 2, 3 
and 4 dimensions, in order to understand these interrelations. Results from simulations of model 
glass formers show that the degree of the breakdown of the SE relation decreases with increasing 
spatial dimensionality. The breakdown itself can be rationalized via the difference between the ac- 
tivation free energies for diffusivity and viscosity (or relaxation times) in the Adam-Gibbs relation 
in three and four dimensions. The behavior in two dimensions also can be understood in terms of a 
generalized Adam-Gibbs relation that is observed in previous work. We calculate various measures 
of heterogeneity of dynamics and find that the degree of the SE breakdown and measures of hetero- 
geneity of dynamics are generally well correlated but with some exceptions. The two dimensional 
systems we study show deviations from the pattern of behavior of the three and four dimensional 
systems both at high and low temperatures. The fragility of the studied liquids is found to increase 
with spatial dimensionality, contrary to the expectation based on the association of fragility with 
heterogeneous dynamics. 



I. INTRODUCTION 



The dramatic slowdown of dynamics upon cooling 
glass-forming liquids towards the glass transition is de- 
scribed by the temperature dependence of a number of 
transport coefficients, and relaxation time scales, which 
include the shear viscosity (77), the translational diffu- 
sion coefficient (D), the rotational correlation time (t c ) 
and the structural relaxation time r a obtained from the 
long time decay of density correlation functions. The 
Stokes-Einstein (SE) relation [IH3] in its original formu- 
lation relates the translational diffusion coefficient (D) 
of a macroscopic, or Brownian, probe particle - a single 
particle property - to the shear viscosity (77) of the liquid 
- a collective property at a temperature T : D = ^^-^ 
where m is the mass and R is the radius of the particle, T 
is the temperature of the liquid and the factor c (= 6 or 
4) depends on the (stick or slip) boundary condition at 
the surface of the Brownian particle [l| . Although a hy- 
drodynamic relation, the SE relation is known to be ap- 
plicable even for the self diffusion of the liquid particles at 
high temperatures Q. However, several experiments and 
simulation studies in the last three decades [5l-[26L [28l - [3l| 
have conclusively shown that at low temperatures in su- 
percooled liquids, the SE relation breaks down, such that 
self diffusion coefficient is much larger than what one may 
infer, using the SE relation, from the value of the viscos- 
ity. Although many works claim that the break down of 
the SE relation occurs around the mode coupling temper- 
ature T c [27| . more recent results (e.g. (29[) suggest that 
the breakdown occurs in the vicinity of a much higher 
temperature, the onset temperature T onset , at which as- 



pects of slow dynamics including stretched exponential 
relaxation (see below) begin to be manifested [32l434| . 
There is a commonly held view that the breakdown of 
the SE relation is a manifestation of dynamical hetero- 
geneity (DH), which is characterized by various indica- 
tors we discuss below, such as the dynamical susceptibil- 
ity (X4)[35|||3, the Kohlrausch- William- Watts (KWW) 
exponent (Pkww) quantifying the stretched exponential 
decay of correlation functions, and the non-Gaussian pa- 
rameter «2 that measures the deviation of particle dis- 
placements from a Gaussian form. We note, however, 
that there is no agreement on the nature and origin of 
heterogeneity among the theories proposed to explain 
the SE breakdown (e.g., dynamical facilitation 41 1, the 
random first order transition theory [4214441 ] . the mode- 
coupling theory [45j . the shear transformation zone the- 
ory [46j] and the obstruction model [13] )• The goal of the 
present work is to interrogate the manner in which the SE 
breakdown changes with the spatial dimensionality, and 
the correlation of SE breakdown with characteristics of 
heterogeneous dynamics, and fragility. Since the spatial 
dimensionality enters explicitly in some of the theoretical 
descriptions mentioned above, this study is expected to 
provide valuable information about the validity of these 
theories. As we discuss below, our investigations reveal 
surprises whose rationalization should lead to a better 
understanding of the key features of glassy dynamics. 

The rest of this paper is organized as follows: In Sec. 
II, we provide a discussion of the SE relation and its 
relation with dynamical heterogeneity and fragility, that 
forms the necessary background for our work. In Sec. 
Ill we provide details regarding the models studied and 
the simulation procedure. Sec. IV contains the results of 
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our study. Sec. V contains a discussion and conclusions 
arising from our work. 



II. STOKES-EINSTEIN BREAKDOWN, 
DYNAMICAL HETEROGENEITY AND 
FRAGILITY 

One of the earliest theories of the breakdown of the SE 
relation is due to Hodgdon and Stillinger who envis- 
age the highly viscous supercooled liquid to be composed 
of sparse "fluid-like" regions of low viscosity in a matrix 
of a more viscous fluid. Thus, both D and n are space- 
dependent. By calculating the viscous drag force on a 
diffusing particle in the "fluid-like" region, they showed 
that the local drag force decreases from the Stokes' value, 
and hence the local diffusion coefficient increases. If one 
uses the bulk viscosity (which is dominated by the more 
viscous regions) in the SE relation, one finds a break- 
down of this relation. In order to apply their model to 
realistic systems and to explain the difference between 
the behavior of translational and rotational diffusion co- 
efficients, Hodgdon and Stillinger had to impose certain 
special properties on their model. However, Tarjus and 
Kivelson [l7[ later argued that the mere existence of do- 
mains with different local properties is a sufficient condi- 
tion for the SE breakdown. They considered a liquid in 
which there are domains (of unspecified nature) of size 
L with a size distribution p(L) such that p(L)L 2 dL is 
the probability of finding a molecule in a domain of size 
between L and L + dL. They assumed that the local re- 
laxation time t(L) in a domain is size dependent and the 
measured a relaxation time (r Q ) is the average of t(L): 



relaxation times. Blackburn et al. [T(| argued that the 
translational diffusion coefficient D and the rotational 
correlation time t c measure different moments of this dis- 
tribution, thus causing the SE breakdown: 

D cx (-), t c = (t) 



Dt c cx (t)(-) 

T 



(r)(-) 



1 for S function distribution, normal SE 
> 1 SE breakdown 



_ (3) 
La Nave et al. showed, [48[ using the potential en- 
ergy landscape framework, that in the 3D Kob- Andersen 
model the product < D >< 1/D > inde ed g rows as 



the temperature decreases. Swallen et al. [12[ argued 
that since Dt c increases with decreasing temperature, 
Eq. [3] implies that the distribution of relaxation time 
should be broader at lower T. 

The stretching exponent Pkww provides a measure of 
the width of the distribution of relaxation times. Let 
p(r) denote the distribution of relaxation times r where 
each of the local relaxation function is exponential with 
relaxation time r. The overall correlation function <f)(t) 
is empirically given by a stretched exponential: 



4>{i) = / drp(r) exp(-i/r) = cx V {-(t/T KWW f KWW ) 



Using simple mathematical identities, it can be shown 
[5fjl ] that the n-th moment (r n ) of the distribution p(r) 
is given by 



t l oc exp(E(L)/k B T) 

r Q = (r L > oc / p(L)cxp(E(L)/k B T)L 2 dL (1) 
Jo 

In their picture, the SE relation is valid inside a do- 
main. However, the translational diffusion involves pas- 
sage through many domains. This is why the average D 
is different from the prediction of the SE relation. As- 
suming that a diffusing particle performs a random walk 
across domains and D(L) changes abruptly at interfaces, 
and neglecting a term involving the gradient of D, one 
obtains: 

D(L) cx T/?](L), n(L) cx cxp(E{L)/k B T) 

D = (D(L)) cx / p(L)D(L)L 2 dL 
Jo 

/>oo 

cx / p{L)cxp(-E(L)/k B T)L 2 dL (2) 
Jo 

Unless p(L) is a (5-function, Eqs. [T] and [2] lead to a viola- 
tion of the normal behaviour, Dt — constant. 

This picture of the heterogeneity may also be inter- 
preted in terms of the existence of a distribution of local 
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Using the above formula, it is easy to show that the rel- 
ative variance, which provides a measure of the width of 
the distribution, is given by: 



(t 2 ) - (t) 2 Pkww^ 
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(5) 



This equation implies that the relative variance depends 
only on (3kww smd increases monotonically as Pkww 
decreases. 

The above analysis suggests the following interpreta- 
tion of the SE breakdown. As the temperature decreases, 
DH (of unspecified nature) develops which leads to the 
existence of a distribution of local relaxation times, the 
width of which increases as T is decreased. This is mani- 
fested as (i) the lowering of Pkww and (ii) the SE break- 
down. Since both the SE breakdown and the lowering of 
Pkww from 1 are manifestations of DH, they should oc- 
cur simultaneously according to this interpretation. 
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We should, however, note the following points that ar- 
gue against this interpretation: (i) the lowering of Pkww 
does not prove the existence of a distribution of relax- 
ation times (Hoj i.e. liquids at low temperatures can be 
dynamically homogeneous with a single relaxation time 
but exhibit an inherently non-exponential decay of corre- 
lation functions and (ii) some experiments on OTP seem 
to suggest that Pkww remains constant in the relevant 
low-temperature range (See [l2[ and references therein). 
In that case the width of the distribution does not change 
with temperature, and Eq. [3] can not explain the ob- 
servation that the degree of the SE breakdown becomes 
progressively larger as the temperature is decreased. 

The DH can be quantified by more direct indicators. In 
the present study, we have used the following quantities 
to characterize the DH. 

a. The dynamical susceptibility X4-' The dynamical 
susceptibility X4,(i), which is the integral of the four-point 
correlation function gi{r,t), measures the fluctuation in 
the two-point correlation function q(t) (overlap function) 
[39[ | . Thus Xi(t) is a direct measure of the DH. For super- 
cooled liquids, Xi{t) shows a peak at a time proportional 
to the a relaxation time J2f|. The peak value of the dy- 
namical susceptibility - X4 - is a direct measure of the 
volume of space correlated during structural relaxation 
j40l | . Experiments and simulations, typically carried out 

in three dimensions, show that X^^iT) grows monoton- 
ically as the temperature T is lowered, which indicates 
that larger regions of space are dynamically correlated 
at lower temperature i.e. the DH is more prominent at 
lower temperature. 

b. The stretching exponent Pkww : The stretching 
exponent Pkww is a measure of how non-exponential 
the decay of a time correlation function is. According to 
the interpretation of the non-exponential decay of the 
density-correlation function as a manifestation of het- 
erogeneous dynamics, a lower value of Pkww implies 
stronger DH. 

c. The non- Gaussian parameter 0:2 (i)." This param- 
eter quantifies the deviation of the distribution of particle 
displacements in time t from the Gaussian form expected 
for spatially homogeneous dynamics. The non-Gaussian 
parameter, calculated from the second and fourth mo- 
ments of the distribution of the displacements, exhibits 
a peak at a characteristic time that increases as the tem- 
perature is decreased. The peak value of ai provides a 
characterization of the spatial heterogeneity of the dy- 
namics. 

d. Fragility: Fragility is a material parameter that 
measures how rapidly the viscosity (or relaxation time) 
of supercooled liquids increases as the temperature de- 
creases. Bohmer et al. j5lj in their extensive compila- 
tion of available data found that the kinetic fragility m 
determined from the slope of the Angell plot at T g has 
negative correlation with the stretching exponent Pkww, 
i.e., smaller Pkww implies a higher fragility. This con- 
clusion remains the same for both isobaric and isochoric 



fragilities [52j . Such a relation has also been obtained 
theoretically within the framework of random first order 
transition theory [43| . Based on this correlation between 
fragility and PkwWi fragility can be considered to be an 
indicator of DH - more fragile systems are more hetero- 
geneous. 

However, there are evidences against this correlation 
as well: (i) In experiments on supercooled water con- 
fined in nano-pores, Chen et al. [22j observed a fragile 
to strong transition and found fractional SE relations in 
both the regimes. However, the breakdown exponent is 
closer to 1 (= 0.74) for fragile water than that for strong 
water (= 0.67). (ii) Similarly, the breakdown exponent 
computed by Jung et al. |4l| using dynamical facilita- 
tion theory was also closer to 1 (= 0.73) for a fragile 
glass-former model than for a strong glass-former model 
(= 0.67). (hi) In [53[ it was shown that the strong cor- 



relation found by Bohmer et al. 51| between the kinetic 
fragility and Pkww becomes much weaker if subgroups 
{e.g. only simple or complex molecul ar g lass-formers) are 
considered, (iv) Dyre also claimed [54|, on the basis of 
experiments on simple, organic glass-forming liquids that 
no clear correlation is present between these two quan- 
tities, (v) In the simulation of ST2 water, Becker et al. 
[23| observed a breakdown of the SE relation at low T in 
both strong and fragile regimes, with the breakdown ex- 
ponent nearly the same for both strong and fragile water, 
(vi) Vasisht and Sastry [55j . similarly find in simulations 
of silicon that the SE breakdown becomes apparent in 
the high temperature, high density liquid, and the break- 
down exponent is the same in both the high temperature 
(fragile) and low temperature (strong) liquids. 

We emphasize, however, that the above picture is 
based on experiments and simulations predominantly in 
three dimensions and does not, a priori, tell us what to 
expect in other spatial dimensions. In the present study, 
we aim to understand the inter-relations among DH, the 
SE breakdown and the fragility in other spatial dimen- 
sions by studying model glass-forming liquids in 2, 3 and 
4 dimensions and by considering both the SE breakdown 
and direct measures (%4, Pkww, a 2) of DH and fragility. 

We end this section by remarking on the study of 
SE relation in two dimensions. Strictly speaking, trans- 
port coefficients such as 77 and D are not well-defined in 
infinitely large two-dimensional systems at equilibrium 
due to the presence of long-time tails [56|, [53] in cor- 
relation functions appearing in the Grccn-Kubo formu- 
lae for these transport coefficients, and for similar rea- 
sons the the use of the Stokes relation is questionable 
[58| . However, these effects are not important for the 
length and time scales considered in our 2D simulations. 
We find well-defined diffusive behavior (mean-square dis- 
placement proportional to time) at long times in all our 
2D simulations, from which the diffusion coefficient D 
can be obtained without any ambiguity. This is consis- 
tent with the results of several existing MD simulations 
of 2D model systems [59j - [62j in which the diffusion co- 
efficient and the viscosity have been computed and the 
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validity of the SE relation has been examined |6C 



III. SIMULATION DETAILS 

In the present study we perform NVT MD simulations 
for the following models : (a) the Kob- Andersen binary 
mixture at the canonical 80 : 20 composition |49| in 4, 3 
and 2 dimensions (denoted by 4D KA, 3D KA, 2D KA 
respectively) at number densities p = 1.60 (4D), p = 
1.20 (3D and 2D); (b) the modified Kob- Andersen model 
(denoted by 2D MKA) at a different composition 65 : 35 
[631 ] at the number density p = 1.20; (c) the 50:50 binary 
mixture of purely repulsive soft spheres with potential 



V{r) 



-10 



in 3 and 2 dimensions (denoted by 3D R10 



and 2D R10 respectively) at a number density p = 0.85 
[HI . The details of the potentials and units are described 
in the corresponding references. The integration time 
step was in the range dt <G [0.001 — 0.006] depending on 
the temperature. An algorithm due to Brown and Clarke 
[65| was used to keep the temperatures constant. System 
sizes were (1) N = 1500 for 4D KA; (2) N = 1000 for 
3D KA; (3) N = 1000 for 2D KA; (4) N = 2000 for 2D 
MKA; (5) N = 10000 for 3D R10 and (6) N = 2048 for 
2D R10. Runlcngths at each temperature were at least 
100r Q (the a relaxation time defined below). 

We have computed the following measures of time 
scales : 

1. Translational diffusion coefficients Da of one type 
of particles measured from the mean squared dis- 
placement (MSD) of that type of particles. 

2. a relaxation times estimated from the time taken 
to decay to 1/e of the initial value of (a) the overlap 
function q(t) [36M391 16(| . (b) the intermediate scat- 
tering function F(k,t), and (c) the self part of the 
intermediate scattering function F s A{k,t) for one 
type of particles. The value of k corresponds to the 
first peak of the partial structure factor of 
one type of particles. 

Our procedure for computing the overlap function 
is defined in (67j . The intermediate scattering func- 
tion F(k, t) and its self part F s (k, t) are defined as 

N N 

F (k,t) = ]T5>-^ (t) -^ (0)) 
i=l j=l 

N 



F s (k,t) = ]T, 



-zk-(fi(t)-fi(0)) 



(6) 



3. The shear viscosity (rj) computed from equilibrium 
simulations using the Einstein and the Green-Kubo 
relations. The Green-Kubo relation for the shear 
viscosity is given by the integral of the auto corre- 
lation function of the stress tensor P a p{t) : 



V 

ksT Jo 



dt < P a p(t)P a p(0) > 



^ t N N N \ 

P a p{t) = y \^2PiaP^/m + ^2^2r ija f ij0 \ (7) 
\l=l 1=1 j>i J 

where ry = \r t - r}| and = - ^0^ and a, ft € 

(x, y, z) denotes Cartesian components. 

The shear viscosity can also be computed from the 
corresponding Einstein relation as the slope of the 
Helfand moment at long time and in linear regime 



1 



'/ 



Vkr>T *->oo 



q < (A aP (t) - A aP (oyf > (8) 



2t 



where A a p(t) is the Helfand moment and related 
to the stress tensor as 



dA aP (t) 
dt 



P a p{t)V 



Since the shear viscosity is a collective property of 
TV particles, numerical accuracy is a big issue when 
computing shear viscosity using either Eq. [7] or 
Eq. [5J To improve numerical accuracy we take 
average over different components. We have found 
that values from the two methods mutually agree 
well and report here the shear viscosity obtained 
from the Einstein method (denoted by rjEinstein)- 



IV. RESULTS 

In this section, wc describe our results for the SE 
breakdown, fragility and different measures of DH in 2, 
3 and 4 dimensions. In Sec. IV (A) we show diffusiv- 
ity, viscosity and relaxation time data that demonstrate 
the SE breakdown. We further show that in 3 and 4 
dimensions, the breakdown exponent can be understood 
in terms of the Adam-Gibbs relation for diffusivity and 
relaxation time. In 2 dimensions, the Adam-Gibbs rela- 
tion is not valid as we have previously shown [68|. We 
show here that nevertheless, the the SE breakdown in 
two dimensional systems can be rationalized in terms of 
a generalized Adam-Gibbs relation observed to hold in 
these systems. In Sec. IV (B) we describe various mea- 
sures by which heterogeneity can be quantified, and the 
correspondence of the degree of heterogeneity with the 
degree of the SE breakdown. In Sec. IV (C) we show that 
fragility increases with increasing spatial dimensionality, 
contrary to expectations that the fragility is correlated 
with the degree of heterogeneity of dynamics. We show 
however that the variation of fragility with dimensional- 
ity can be understood in thermodynamic terms, through 
the evaluation of the configurational entropy and proper- 
ties of the distribution of local energy minima or inherent 
structures I69II. 
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The dimension dependence of the SE 
breakdown 



The breakdown of the SE relation in three dimensions 
in the Kob-Andersen model has been reported earlier 
[TBI [29l [30| . Our data agree reasonably well with those 
of previous studies. Here we study the dependence of 
the SE relation on the spatial dimension by considering 
the following models which are defined in the previous 
section: (i) 2DKA, (ii) 2DMKA (iii) 2DR10, (iv) 3DKA 
and (v) 4DKA. According to the SE relation, the quan- 
tity Drj/T should be independent of the temperature T 
- deviations of this ratio from a constant value as the 
temperature is changed would constitute a violation of 
the SE relation. Since the viscosity i] is difficult to cal- 
culate in a molecular dynamics simulation, the tempera- 
ture dependence of Dr a or Dr a /T has been considered 
in most existing numerical studies to investigate the va- 
lidity of the SE relation. These choices correspond to a 
replacement of rj by r Q assuming either that r Q oc r\ or 
T a oc rj/T. The justification for these substitutions ob- 
tain from either the Maxwell relation (77 = Goqt) or the 
relationship between relaxation times and the diffusion 
coefficient in the diffusive regime (t -1 = Dq 2 = mkB ^ X. ■ 
this argument assumes the validity of the SE relation). 
From various experimental and simulation results, [70l — 
[72j|, T a oc r//T. Consistently with previous observations, 
we find that r Q oc rj/T provides a very good description 
of the low-temperature data over a fairly large range of 
rj and r Q , as shown in Figfl] using viscosities calculated 
from the Green-Kubo formula and the Einstein relation 
for the 3DKA model. We therefore use the following two 
representations to analyze the degree of SE violation from 
the data. 

i In a Da vs. t u or ^ plot (Da is the calculated diffusion 
coefficient of particles of type A), by fitting the data to 

a power law of the form Da oc t~^, or Da oc (^) ^ we 
test if the effective exponent £ is different from unity, 
the prediction of the SE relation. 

ii From the T dependence of DAT a or D ^ , we test if 
these quantities are constant (SE relation obeyed) or 
become T dependent (breakdown of the SE relation). 

The breakdown of the SE relation in 2, 3 and 4 di- 
mensions is shown in Figs. [2] (2DKA), [3] (2DMKA), H 
(2DR10), (S^3DKA), (B{3DR10) and (4DKA). In most 
of the cases (excepting 2DKA and 3DR10) we show data 
ranging from high temperatures, well above the onset 
temperature T onset , to well below T onset . A change from 
Arrhenius to super-Arrhenius temperature dependence 
of relaxation times defines the onset temperature of slow 
dynamics [32T - I341 ] T onset . The onset temperature is close 
in all cases to the temperature at which one observes the 
breakdown of the SE relation, Tseb as described below. 
All the simulations reported here are performed at tem- 
peratures above the mode coupling temperatures T c for 
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3D KA, N=1000,Density=1.20 
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FIG. 1: a relaxation times r obtained from (i) the overlap 
function (r qt , (ii) the self part (rp s Akt) of the intermediate 
scattering function F s A(k* ,t) of one type (A) of particles 
and (iii) from the full (rFkt) intermediate scattering func- 
tion F(k*,t) plotted against viscosity showing that r oc rj/T 
is a good description of data at low T in the 3D KA model. 
Systematic deviations are seen at high T. k* is at the first 
peak of the partial structure facture SAA(k). 



the respective models, and the T$eb we obtain are well 
above the T c values. 

one below the other in single column instead of messing 
with two column figures. 

From the Da vs. r a or ^ plots, we see that in all 
spatial dimensions, the low T data follow a fractional 
power-law relation, indicating a breakdown of the SE re- 
lation. For the same model (KA) in different dimen- 
sions the power-law exponent £ is closer to 1 in higher 
dimensions, indicating that the SE breakdown is weaker 
in higher dimensions |73l . [74| ■ The difference in the expo- 
nent in two and three dimensions however is very small 
and essentially negligible. Similarly for the RIO model, 
the exponent is marginally, but negligibly, higher in three 
dimensions. The MKA model in two dimensions has an 
exponent that is very close to the KA model in the same 
dimension. All the model results taken together, the 
breakdown exponent is different for different models in 
the same dimension. 

Further, as the temperature increases, there is a clear 
change of the exponent value in Da vs. t q or ^, plots 
indicating a qualitative difference between high T and low 
T behaviours. Surprisingly, in two dimensions, in both 
the models for which high T data are shown (2DMKA 
and 2DR10), the high T exponent from Da vs. r a plots is 
higher than 1. However, in 3 dimensions, the SE relation 
Da oc ^ is recovered at high T as expected. Finally, as 
we go to still higher dimension (D = 4), the Da vs. r a 
plot shows the expected relation (Da oc r" 1 ) at high T. 

The above correlation between Da and r a is also re- 
flected in the T dependence of D AT a . Since the dif- 
fusion coefficient decreases but the relaxation time in- 
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FIG. 2: Plot showing the breakdown of the SE relation in 
the 2D KA model at low temperatures. Top: Da vs. t plot 
(Here r is the a relaxation time from the overlap function 
q(t)). Bottom: T dependence of Da t. The data shown follow 
a fractional SE relation. 



creases as the temperature decreases, there is a competi- 
tion between two opposing effects. The respective rates 
of increase and decrease with T exactly cancel each other 
only if the exponent is 1. Since in 2D, the exponent is 
never 1, the quantity Dat u goes through a minimum and 
approaches a constant value only at very high tempera- 
tures. In 3D and 4D, Dat becomes constant at high 
T as expected and it increases as T decreases at low T, 
indicating a breakdown of the SE relation. 

We have also considered the behaviour of Da/T vs. 
r a and the T dependence of D y° . This corresponds to 
assuming that r Q cx r\. The Da/T vs. r a plots show a 
fractional SE relation at low T in all dimensions. The 
values of the breakdown exponents from the Da/T vs. r a 
plot are slightly different from those obtained from the 
Da vs. r a plots. However, they show the same trend: 
the magnitude of the exponent is closer to 1 at higher 
dimensions, indicating that the SE breakdown is weaker 
at higher dimensions. The estimates of the breakdown 
exponents arc summarized in Table HI 



TABLE I: Estimates of the magnitude of the SE breakdown 
exponents in different spatial dimensions D. Notations: (a) 
^se _ breakdown exponent obtained from Da vs. t q or 
Y plots; (b) £ AG = ratio of slopes from AG plots using Da 
and r a ; High T exponents are obtained from Da vs. r a plots. 







Low T Exponents 




D 


Model 


^SE 




High T exponents 


2 


2DR10 


0.75 




1.18 


2 


2DKA 


0.82 






2 


2DMKA 


0.84 




1.50 


3 


3DR10 


0.752 






3 


3DKA 


0.83 


0.85 


« 1 


4 


4DKA 
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Figs. [21 HI [S] and [7] also show that the temperature of 
SE breakdown (Tseb), estimated as the point of inter- 
action of the low T and the high T fits in Da vs. r a 
or ^ plots, is close to the Arrhenius to non-Arrhenius 
cross-over temperature T onset in all dimensions. Thus, 
the SE breakdown occurs at a temperature that is much 
higher than the divergence temperature T c of the mode- 
coupling theory estimated from a power-law fit of the 
T-dependence of the relaxation time. 

The fractional SE behaviour at low T can be ratio- 
nalized by considering the different dependence of the 
diffusion coefficient and the a relaxation time on the con- 
figuration entropy (S c ). The Adam-Gibbs (AG) relation 
X = Xq exp(^^) , if it is valid, provides a way to test 

this hypothesis quantitatively (X is r a or (Da)' 1 in the 
present study). Fig. [8] shows that the AG relation is 
valid in the 3DKA model (top row) and the 4DKA model 
(bottom row). We see that the slope of the (Da)' 1 vs. 
(TiSc) -1 plot is different from that of the r Q vs. (TSc)' 1 
plot. Table U shows that the observed fractional SE ex- 
ponent at low T can be interpreted as the ratio of the 
slopes in the AG plots in Fig. [5J 

In the two dimensional models one sees a deviation 
from the AG relation at low temperatures. But the be- 
haviour can still be rationalized by considering the gen- 
eralized Adam-Gibbs relation that is a good description 
of the data [68|]. Thus, we can write 



In Da 
InT 
In D A 



In D 

In t 4 

In D 

otp 



Ad_ 

TS C/ 

TSc) 

Ad ' 



ln^ 



(9) 



Eq. El provides the relationship between D and r 
based on each obeying a generalized Adam-Gibbs rela- 



7 




FIG. 3: Plots showing the breakdown of the SE relation in the 2DMKA model. Left: Da vs. t plot (Here r is the a relaxation 
time from the overlap function q(t)). Right: T dependence of Dat. The low T data follow a fractional SE relation. A clear 
change of exponent occurs at high T in Da vs. t plot, although the high T exponent is bigger than 1. The change of slope 
occurs at a temperature Tseb which is close to T onae f {Tseb estimated as the point of intersection of high T and low T 
fits; T onse t is the onset temperature of slow dynamics.) Also shown is the dependence of D on r according to the generalized 
Adam-Gibbs relation discussed in the text. 
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FIG. 4: Plots showing the breakdown of the SE relation in the 2DR10 model. Left: Da vs. t plot (Here r is the a relaxation 
time from the overlap function q(t)). Right: T dependence of Dat. The low T data follow a fractional SE relation. A clear 
change of exponent occurs at high T in Da vs. t plot, although the high T exponent is bigger than 1. The change of slope 
occurs at a temperature Tseb which is close to T onse t. (Tseb estimated as the point of intersection of high T and low T 
fits; Tonset is the onset temperature of slow dynamics.) Also shown is the dependence of D on r according to the generalized 
Adam-Gibbs relation discussed in the text. 



tion, which arc shows for 2DMKA and 2DR10 models in 
Figs. El |H It is seen that indeed the low temperature 
data as well as some part of the high temperature data 
are well described by this form. It must be noted that 
the Eq. [HI is not a fractional SE relation, and thus, the 
description of the data for these two dimensional systems 
through a breakdown exponent is an approximation. 

To summarize, results presented here show that: (i) 
the SE breakdown in weaker in four dimensions than in 
three dimensions which is consistent with earlier works 
[73l [74| . (ii) The breakdown exponent can be rational- 
ized from the different scaling of the diffusion coefficient 
and the relaxation time with the configurational entropy, 



either via the AG relation (three and four dimensions) 
or a generalized AG relation (two dimensions) . (iii) The 
behaviour in two dimensions is more complicated, dis- 
playing no temperature regime where the SE relation is 
valid, consistently with earlier work (75|. (iv) The SE 
breakdown temperature Tseb is, in all cases, very close 
to the onset temperature T onset . 
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FIG. 5: Plots showing the breakdown of the SE relation in the 3DKA model. The a relaxation times r are computed from 
the (i) overlap function (r q t) (ii) F(k,t) at the peak of S(k) (rpkt) These two measures are mutually proportional and are 
used interchangeably. Left: Da vs. t and ^ plot. r q t, Tpkt(k*) are multiplied by constant factors to match all data sets at 
low temperature. Right: T dependence of Dat and 3-^-. The low T data follow a fractional SE relation. A clear change 
of exponent occurs at high T in the Da vs. t plot. The high T exponent (=-1) expected from the SE relation is obtained 
from the Da vs. ^ plot. The SE breakdown occurs at a temperature Tseb closer to the onset temperature T onse t than T c 
(Tseb estimated as the point of intersection of high T and low T fits; T onse t is the Arrhenius to non-Arrhenius cross-over 
temperature; T c is the mode coupling transition temperature. All data points shown here are at T > T c .) 




FIG. 6: Plot showing the breakdown of the SE relation in the 3D RIO model at low temperatures. Left: Da vs. t plot (Here 
r is the a relaxation time from the overlap function q(t)). Right: T dependence of Dat. The data shown follow a fractional 
SE relation. 
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FIG. 7: Plots showing the breakdown of the SE relation in the 4DKA model. Left: Da vs. t plot (r here is the a relaxation 
time from the overlap function q(t)). Right: T dependence of Dat. The low T data follow a fractional SE relation. A clear 
change of exponent occurs at high T in Da vs. t plot. The high T exponent (=-1) expected for a homogeneous (Gaussian 
distribution of particle displacements) system is obtained from Da vs. t plot. The SE breakdown occurs at a temperature 
Tseb closer to the onset temperature T „ se t than T c . (Tseb estimated as the point of intersection of high T and low T fits; 
T onse t is the Arrhenius to non-Arrhenius cross-over temperature.) 



3D KA, N=1000, Density=1.20 
Ratio of slopes = 0.85 













liD 


' = 1.49/TS c +2.64 





t 






1m 


=1.75/15, -10 



3.0 4.0 

crsr' 



4DKA,N=1500,Density=1.60 
4_ Ratio of slopes = 0.9 



. 10 : 



M) A =W/TS [ t2.97 
]'nl=3Jl/TS -0.72 



iTS.j" 



FIG. 8: Left: The Adam Gibbs (AG) relation in the 3DKA model using as the dynamical quantities Da, t. Here r is the a 
relaxation time from the overlap function q(t). Right: The AG relation in the 4DKA model using as the dynamical quantities 
Da, r. The slopes are different for Da and r indicating that the diffusion coefficient has a different dependence on the 
configuration entropy than the a relaxation time. The fractional SE exponent at low T can be interpreted as the ratio of the 
slopes (Table H]). 




FIG. 9: The time dependences of the dynamical susceptibility x-i(t) an d of the non-Gaussian parameter ce2,A(t) in the 4DKA 
model. Peak values ~x^ a and o£f\ " extracted from such time dependence are described further below. 
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FIG. 10: ieft: T dependences of the peak height x? 6 " °f the dynamical susceptibility %4(t) in the 2D RIO model. .Right: T 

peak 
V 2,A 



dependences of the peak height af e jj fe of the non-Gaussian parameter 0:2, a(£) in the 2D RIO model. 
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FIG. 11: £eft : T dependences of the peak height y^ ak of the dynamical susceptibility Xi{t) in the 2DMKA model. .Right : T 



dependences of the peak height af 6 ^* of the non-Gaussian parameter 02, a (t) in the 2DMKA model. 





FIG. 12: ieft : T dependences of the peak height y^ ah of the dynamical susceptibility Xi(^) m t ne 3DKA model. Right : T 
dependences of the peak height c£^°^ of the non-Gaussian parameter ai,A{t) in the 3DKA model. 
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FIG. 13: ieft : T dependences of the peak height y^ ah of the dynamical susceptibility Xi(t) m t ne 4DKA model. .Right : T 
dependences of the peak height af e jj fe of the non-Gaussian parameter 02, a (t) in the 4DKA model. 




FIG. 14: T dependence of F sA kt in 2DR10, 2DMKA, 3DKA and 4DKA models. Also shown are fit curves to the stretched 
exponential form. Horizontal lines indicate the range of y values within which the fitting is performed. 
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Dynamical heterogeneity and the breakdown of 
the Stokes-Einstein relation 



Here we consider three measures of dynamical het- 
erogeneity, namely the dynamical susceptibility Xii the 
non-Gaussian parameter 012 and the stretching exponent 
Pkww- We examine to what extent these measures of 
heterogeneity correlate with the degree of the breakdown 
of the SE relation. 

Both X4 and ct2 have been extensively studied as in- 
dicators of dynamical heterogeneity (e. g. [37i - [39l |7(|. 
These quantities are functions of time, and exhibit a 
maximum value at a characteristic time that is close to 
r a [39| in the case of X4 and to t* ~ (D/T)^ 1 in the 
case of «2 uM- The time dependence of X4 and 012 are 
shown in Fig. |H1 Peak values y^ ak and a_2 peafe extracted 
from such time dependence are described further below 
as measures of heterogeneity. 

The temperature dependence of the peak of the dy- 
namical susceptibility X4 6 ° and that of the non-Gaussian 
parameter oP^ ak in different models are shown in Figs. 
[POl [TTl [T2l and rrjl Clearly, in all dimensions, the peak 
heights X4 £ak and Oi^ eak grows as the temperature de- 
creases thus indicating that the dynamics is heteroge- 
neous in the lower part of the temperature range that we 
study. 

As a third direct measure of the degree of dynami- 
cal heterogeneity, we estimate the stretching exponent 
Pkww- We have computed Pkww from the self inter- 
mediate scattering function F S Akt for k at the peak of 
the structure factor. Fig. 1141 shows the temperature de- 
pendence of F s Akt{k,t) for different models and the fits 
to stretched exponential form from which the stretching 
exponents Pkww are obtained [7?} ■ 

In order to compare the degree of the SE breakdown 
with that of dynamical heterogeneity, we plot X4 6afc ; 
oUr^A* and /3kww f° r the different models against the 
relaxation time r Q . We expect that the system show- 
ing a stronger SE breakdown at a given temperature (or 
relaxation time) will also show a larger degree of hetero- 
geneity. The data in Fig. [T5]for \ v ^ ak shows that indeed, 
the four dimensional system by and large has smaller val- 
ues of Xi eak than the three dimensional system, although 
the values are close especially at low temperatures. This 
is consistent with the larger SE breakdown observed in 
three dimensions. Two dimensional systems show sub- 
stantially larger values of x v ^ ak than three dimensional 
systems for most of the studies points. While this is quali- 
tatively consistent in principle with the marginal increase 
in the degree of SE breakdown in two dimensions, quan- 
titatively the difference in xl*^ is large and surprising. 



In similar fashion, oSf^ values show an increase in 
going from four to three dimensions. However, two di- 
mensional systems show smaller values of oE^"^ than ei- 
ther three or four dimensional systems. This is contrary 
to the relative degrees of SE breakdown in these differ- 



ent systems. Such an inconsistency may reveal either a 
shortcoming of o£f£* as a measure of heterogeneity, of 
another aspect of the peculiar behavior of two dimen- 
sional systems. Further investigations are necessary to 
understand this behavior. 

Finally, we find that the stretching exponents Pkww 
depend systematically on the dimensionality, even 
though the estiamtes are noisy. The Pkww values for 
the four dimensional system (4DKA) lie above those of 
the three dimensional system (3DKA), which in turn are 
larger than those the values for the two dimensional sys- 
tems (2DMKA and 2DR10). Though not strictly consis- 
tent with the degree of SE breakdown as quantified by 
the breakdown exponent, the behavior of Pkww does not 
throw up any surprises. A more detailed study, consider- 
ing also length scales over which dynamics is correlated 
spatially, perhaps along the lines of (7l| . 
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FIG. 15: Comparison of the degree of heterogeneity in differ- 



ent dimensions using the y e ^ a as a measure. Vertical lines 
correspond to Tseb 



Dependence of the fragility on spatial 
dimensions 



As discussed in the introduction, the fragility of a glass 
former has been argued to be correlated with the hetero- 
geneity of a glass forming liquid, and by extension the 
degree of the breakdown of the SE relation. To evalu- 
ate such an expectation, we calculate the fragility of the 
model glass formers studied. We estimate the kinetic 
fragility from the Vogcl-Fulchcr-Tammann (VFT) fit to 
the T-dcpcndence of the a relaxation times: 



r(T) = r exp 



K 



VFT 



Tvi 



1) 



(10) 



Fig. HU shows the T dependence of the a relaxation 
times in the five models we have studied, in a "fragility 
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FIG. 16: Comparison of the degree of heterogeneity in differ- 
ent dimensions using aS cafe 



as a measure. 




FIG. 17: Comparison of the degree of heterogeneity in differ- 
ent dimensions using the Pkww as a measure. 
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FIG. 18: Top: Fragility plot of relaxation times for five mod- 
els in 2,3,4 spatial dimensions. VFT fits to relaxation times 
shown are used to obtain the kinetic fragility Kvft (see dis- 
cussion). The fragility plot employs a "simulation glass tran- 
sition temperature" T g defined as r(T g ) — 10° (reduced unit) 
to scale temperatures. The Kvft values are listed in Ta- 
ble [TTJ Bottom: Fragility plot for four models in 2,3,4 spatial 
dimensions using the inverse of the diffusion coefficient. The 
plots show that systems at higher dimensions are more fragile. 



plot" which shows relaxation times in an Arrhenius plot 
with the temperature scaled by the glass transition tem- 
perature. We use a "simulation glass transition temper- 
ature" T g dchned from r(T g ) = 10 5 (reduced units) in 
order to make such a plot. The bottom panel shows a 
fragility plot for the inverse diffusion coefficients, with T g 
defined from D^ 1 = 10 5 (reduced units). The plots and 
the VFT fits to the relaxation times clearly show that 
the fragility of the liquids is larger at higher dimensions. 
This trend is opposite to the expectation based on that 
of the breakdown of the SE relation, and an assump- 
tion that the two are correlated. Nevertheless, there are 
no strong reasons to expect that fragility is correlated 
with heterogeneous dynamics, in particular with a varia- 
tion of dimensionality. For example, the analysis by Xia 



and Wolynes [43j is restricted to three dimensions, and it 
would be interesting to carry it out for other dimensions 
as well. 

On the other hand, a rationalization of the fragilities 
themselves of the different systems studied can be at- 
tempted using the energy landscape approach in [69| . 
based on the validity of the Adam-Gibbs relation. In 
order for the VFT relation and the Adam-Gibbs relation 
to both hold, one needs the T dependence of TS C to be 
linear: 

TS C = K T - l\ , (11) 
Eq. [TT] defines the thermodynamic fragility Kt- The 
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FIG. 19: T dependence of TS C of five models in 2,3,4 spatial 
dimensions plotted as TS C vs. T/Tk so that the slope is an 
estimate of the thermodynamic fragility Kt (listed in Table 
HI]) . The plot shows that the thermodynamic fragility increases 
with increasing spatial dimensionality. 
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The origin of Kt can be understood from the proper- 
ties of the minima of the potential energy landscape [69J , 
known as "inherent structures" (IS). The configurational 
entropy is the entropy associated with the multiplicity 
of energy minima, and if the density of states of such 
minima with respect to their energy eis is described by 



S c (ei S )/Nk B = a 



-is) 



(13) 



such that a quantifies the total number of inherent struc- 
tures, and a 2 the variance of the distribution, the ther- 
modynamic fragility (under the simplifying assumption 
of the vibrational entropies of all the inherent structures 
being the same) is given by 



K T = Vacr/2. 



(14) 



In other words, Kt is proportional to the total spread 
of the density of states. Fig. ITOl shows the T dependence 
of the configurational entropy for five models. We see 
that the thermodynamic fragility increases significantly 
as the spatial dimension increases. Further, in Fig. I2"0"lwc 
show that the quantity a^fa estimating the width of the 
configurational entropy density increases with increasing 
spatial dimension, which partially explains the increase 
of the thermodynamic fragility with increasing spatial 
dimensions. 

Eq. [12] formally resolves the contributions of configu- 
rational entropy (Kt) and the energy barrier (A) to the 
kinetic fragility. Figs. [19] and [20] shows that the config- 
urational entropy contribution increases with increasing 
dimension thus by itself explaining the increase in the 
kinetic fragility at higher dimensions. The characteristic 
parameters related to fragility for the five models studied 
here are summarized in Table HI] 



SUMMARY AND CONCLUSIONS 



FIG. 20: The configurational entropy density S c (eis) of the 
inherent structures plotted vs. their energy shifted by its 
minimum possible value eis-min defined from S c (eis m i n ) = 



Sc(e IS ) 



0. The lines are fits to the parabolic form: 

- — =yJ^- where e° IS — eismin + o^fa is the IS energy at 
The distribution is broader for 



the peak of the distribution, 
higher dimensions which partially explains the increase of the 
thermodynamic fragility with increasing spatial dimension. 



thermodynamic fragility Kt can be related to the ki- 
netic fragility Kvft by combining the AG relation r = 
r exp(^-), the linear T dependence of TS C (Eq. [11]) 
and the VFT form (Eq. [TO]) to obtain [g3|, 



K 



VFT 



K T /A 



(12) 



In this work we have attempted to analyze the break- 
down of the Stokes-Einstcin relation in systems of spatial 
dimensionality 2, 3 and 4, and its relation to the ther- 
modynamics of the system in the form of the variation of 
configurational entropy, and various aspects of dynamical 
heterogeneity. We summarize here the salient findings of 
this work: (i) We find that in systems in all spatial di- 
mensions, the low temperature relationship between dif- 
fusion and structural relaxation times is well described 
by a fractional Stokcs-Einstein relationship, (ii) The 
high temperature behavior is consistent with the Stokcs- 
Einstein relationship in four and three dimensions, but 
we find that in the two dimensional systems we study, 
the SE relation is not valid at high temperatures also. In 
particular, representing the behavior as a fractional SE 
relationship, we find the exponent to be bigger than 1. 
(iii) We find that the observed breakdown at low tem- 
peratures can be rationalized and understood in terms of 
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the different activation free energies for diffusion and re- 
laxation times in three and four dimensions, and in terms 
of a generalized Adam-Gibbs relation in two dimensions, 
(iv) We find that the crossover from high temperature 
to low temperature behavior (which constitutes a break- 
down of the SE relation for three and four dimensional 
systems) occurs close to the onset temperature of slow 
dynamics, rather than close to the mode coupling tem- 
perature as previously discussed, (v) We find that the 
exponent characterizing the SE breakdown is different 
for different systems even for the same spatial dimen- 
sions (as seen in two and three dimensions) calling into 
question theories which claim universal exponents for a 
given spatial dimension, (vi) We find that the degree of 
breakdown of the SE relation correlates reasonably well 



with measures of dynamical heterogeneity, but with no- 
table exceptions such as the temperature dependence of 
the non-Gaussian parameter ct2- (vii) We find that the 
dependence of fragility on spatial dimensionality is the 
opposite of the degree of the breakdown of the SE rela- 
tion, with fragility becoming larger in higher spatial di- 
mensions. Although surprising, a strong case for a corre- 
lation between fragility and heterogeneity is not a given, 
and arguments for such a correlation in three dimensions 
need to be examined in light of our results for their ap- 
plicability to other dimensions. It would be interesting 
to pursue some of the questions raised by our work in 
higher spatial dimensions as well, such as recently done 
in attempting to identify an upper critical dimension [78| , 
which will hopefully clarify some of the open questions. 
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TABLE II: Fragility-related parameters for the models in different dimensions. Tk is the Kauzmann temperature obtained from 
extrapolating TS C to zero. Kvft is the kinetic fragility obtained from VFT fits to relaxation times. Kt is the thermodynamic 
fragility obtained from the T-dependence of TSc. A is the Adam-Gibbs coefficient. Kag = Kt/A is the kinetic fragility 
estimated from the AG relation. 



Model 


Density 


T K 


TypT 


T 


Kvft 


K T 


A 


K AG = 

K T /A 


4D KA 


1.60 


0.525 


0.530 


0.676 


0.30 


0.972 


3.382 


0.29 


3D KA 


1.20 


0.28 


0.295 


0.402 


0.21 


0.314 


1.79 


0.17 


2D KA 


1.20 


0.477 


0.501 


0.852 


0.11 


0.260 






2D MKA 


1.20 


0.251 


0.214 


0.361 


0.12 


0.166 






2D RIO 


0.85 


0.181 


0.326 


0.453 


0.22 


0.177 







